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. The critical and multicritical behavior of the simple cubic Ising model with nearest-neighbor, 

' next-nearest-neighbor and plaquette interactions is studied using the cube and star-cube approxi- 

mations of the cluster variation method and the recently proposed cluster variation-Pade approx- 
imant method. Particular attention is paid to the line of critical end points of the ferromagnetic- 
paramagnetic phase transition: its (multi) critical exponents are calculated, and their values suggest 
that the transition belongs to a novel universality class. A rough estimate of the crossover exponent 
is also given. 
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C/j ' The Ising model was introduced in the 20's Q for 
describing magnetic transitions and is still the subject 
of a very intense research activity. In particular, the 
Ising model with nearest-neighbor (n.n.), next-nearest- 
neighbor (n.n.n.) and plaquette interactions on the sim- 
^ , pie cubic lattice has been considered as a simple model 
Q for the statistical mechanics of random surfaces [||^ , mi- 
^ O^ , croemulsions ||^, and also as a discretized string action 
(the so-called gonihedric model) |l6|-|l0|. The critical be- 
T-H , havior of the square lattice version of this model has been 
studied for many years and is now well-established 
but very little has been done in the three-dimensional 
case. Previous mean-field calculations H,|l^ have shown 
, that the model exhibits a very rich phase diagram, with 
T— I ' lamellar and ordered bicontinuous |l^ phases, disordered 
^^O , structured and non-structured regions |]l^ , |l5| , and co- 
existing ferromagnetic and paramagnetic phases with a 
first-order wetting transition [|6|. 

The purpose of the present paper is to investigate 
, the critical properties of the model beyond the mean 
' field level ||l^, focusing mainly on the multicritical and 
Ch , crossover properties of the line of critical end points of 
^ ■ the ferromagnetic-paramagnetic phase transition. This 
will be done by means of the cluster variation method 



(CVM) 1 18 1^ in its cube and star-cube Ig^l approxima- 
, tions, and the recently proposed cluster variation-Fade 
■ approximant method (CVPAM) [2l|-|2^. Before turning 
5^ '< to the description of our results, we now give a short 
account of both the model and the method. 

The model is defined by the reduced hamiltonian Ti = 
-f3H 

H = Ji SiSj + J2 ^ SiSj + J3 ^ SiSjSkSi, 

<ij> «ij>> [i,j,k,l] 

(1) 
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where ct; = ±1 is the Ising variable associated to the site i 
of our simple cubic lattice, and sums are respectively over 
the n.n. pairs, n.n.n. pairs and plaquettes of the lattice. 

In terms of the Peierls surfaces separating domains of 
spins with different sign, not only the area but also the 
bending and the intersections of the Peierls interfaces are 
weighted by the couphngs Ji, J2, @]- If (3a, (3c, (3i are 
respectively the energy cost for a plaquette, a bending 
between two adjacent plaquettes and an intersection of 
four plaquettes sharing a common dual link of the Peierls 
interfaces, the relation with the couplings Ji, J2, J3 is 
(3a = 2Ji + 8J2, (3c = 2(J3 - J2), /?/ = -4(J2 + J3) j|. 
The model (|l|) can be considered as a discrete realization 
of a random surface model with an extrinsic curvature 
energy term Recently, the special case in which 

only bendings and intersections are taken into account 
and the area is not weighted at all (that is Ji/(3 = 2k, 
J2//3 = ~k/2, J3/(3 = (1 - k)/2, where /t = is the 
case of phantom surfaces ]2^ ], while k —> 00 represents 
the limit of complete self-avoidance) has been put in con- 
nection with a discretized string model (the so-called go- 
nihedric model) |^ and studied by Savvidy and Wegner 
0. This choice of the couphngs corresponds to a zero 
temperature high degeneracy point where all possible se- 
quences of and "— " planes have the same energy. A 
phase transition has been found in this restricted parame- 
ter space with exponents different from the usual 3d Ising 
exponents [p[-pO[. It has to be observed that this special 
case corresponds to the disorder line J2 = — Ji/4 as 
calculated in the mean-field approximation and that 
in the two dimensional case there is no transition neither 
on the disorder line nor on the line J2 = — Ji/4 [ p9[ . 

We will concentrate on the study of the phase diagram 
around the region where the lamellar, the ferromagnetic 
and the paramagnetic phases coexist, which is close to 
the line J2 — — Ji/4, and we will explain the origin of 
the transition found in . There is a particular phys- 
ical interest in this region due to the extremely low val- 
ues of the surface tension between coexisting phases , 
which is an important property for applications in real 
surfactant systems 

Let us now briefly discuss the methods we are going 
to use. The cluster variation method is a powerful gen- 
eralized mean field theory introduced by Kikuchi fl^ ] 
and then reformulated in a very elegant way [is[| as a 
truncated cluster (cumulant) expansion of the variational 
principle of statistical mechanics. In the cube approxi- 
mation of the CVM one has to minimize the free energy 
density functional 
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where is the contribution of a single cube to the hamil- 
tonian (when sphtting the total hamiltonian H into single 
cube contributions one has to keep in mind that nearest 
neighbour interactions are shared by four cubes and then 
will get a coefficient 1/4 in Hg, , and similarly next-nearest 
neighbour and plaquette interactions will get a coefficient 
1/2), C{x) — x\nx, Pa with a = 8 (4, 2, 1) denotes the 
cube (respectively plaquette, edge, site) density matrix, 
and the sums in the entropy part are over all plaquettes 
(edges, sites) of a single cube (notice that we have not 
assumed any a priori symmetry property for our density 
matrices, and that the plaquette, edge and site matrices 
can be thought of as partial traces of the cube matrix) . 
In the following we shall also use, for the ferromagnetic 
phase only, the star-cube approximation, introduced in 
pO| , where it is described in great detail for the n.n. case 
(inclusion of n.n.n. and plaquette interactions is indeed 
straightforward). For both approximations, the (numer- 
ical) minimization task is greatly simplified by the so- 
called natural iteration method [ p0|j3^ . 

Being an approximate variational theory, the CVM 
yields necessarily classical values of the critical expo- 
nents. In order to overcome this major drawback, one 
of us has proposed the cluster variation-Pade approxi- 
mant method |2|-|2|. The basic idea of the CVPAM is 
that, since the CVM with 7-8 point or larger clusters is 
very accurate at high and low enough temperatures, one 
can extrapolate the results at such temperatures using 
Dlog Fade approximants (see for an application 
of the more sophisticated Adler's methods) in order to 
extract accurate information about the critical behavior, 
i.e. improved critical temperatures and non-classical, pre- 
cise critical exponents. In test applications, indeed, the 
CVFAM has produced results of quality almost compa- 
rable to state of the art Monte Carlo simulations, but 
with a much smaller numerical effort. 

We can describe our results. The ground states of the 
model have been thoroughly investigated in here we 
shall consider the restricted parameter space Ji , J2 and 
K > with J3/J1 = (1 — k)/(4k). In this region, the line 
J2/J1 = —1/4 is always the boundary between ferromag- 
netic and lamellar ground states. 

In Fig. |l| the phase diagram of the model (|l]) as given 
by the CVM cube approximation is depicted in the plane 
J3 — 0. The line separating the paramagnetic and the 
ferromagnetic phase (dashed line) is a second order line, 
with the usual critical Ising transition at Ji = 0.218 
(shifted to Ji = 0.222 by application of the CVPAM 
[21|), to be compared with the best estimate Ji — 0.22165 
[32|. The lamellar phase - here consisting of alternate 
planes of different sign - is separated from the paramag- 
netic and the ferromagnetic phase by a coexistence line 
(solid line) |33| , which is asymptotically close to the line 
J2 — — Ji/4 at low temperature. The second order line 
ends onto the first order one with a critical end point at 
J™d ^ 0.865 ± 0.005 and Jf"'^ = -0.2176 ± 0.0006. In 
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the case of the gonihedric model [|6| with k = 1 , our value 
for the inverse critical temperature is (3c = 0.427 to be 
compared with the value Pc = 0.44 found by Monte Carlo 
simulations jsj and the upper bound f3c = 1-49 calculated 
in 1^. Notice that, since the lamellar- ferromagnetic co- 
existence line is slightly bent toward the lamellar phase 
(and this feature persists at J3 7^ 0), the critical point 
of the gonihedric model is extremely close to our critical 
end point, and this has important consequences on the 
meaning of the exponents that one can define, as we shall 
see. 

The topology of the phase diagram at varying k re- 
mains the same as at J3 = 0, but in the range < k < k*, 
with K* ~ 0.8, there is a tricritical point (see Fig. H) on 
the ferromagnetic-paramagnetic transition line at nega- 
tive values of J2 . The transition line becomes of first 
order before reaching the lamellar phase. 

The critical behaviour suggested by the phase diagram 
of Fig. |l|can be now discussed. The first important point 
is that the critical end point we have found (actually 
a line of critical end points, since n can vary), which 
does not exist in two dimensions |Tl| , must be described 
by critical exponents which differ from the usual three- 
dimensional Ising ones, which apply to the ferromagnetic- 
paramagnetic critical surface. As a consequence, in the 
vicinity of the line of critical end points, which bounds 
the critical surface, it is natural to expect some crossover 
phenomenon. 

In Rcf. we have already calculated with the CV- 
PAM the order parameter critical exponent f3 of the 
gonihedric model (using CVM results up to a temper- 
ature which was less than half the transition temper- 
ature), finding /3 = 0.062 ± 0.003 (together with the 
improved estimate for the inverse critical temperature 
Pc = 0.434), which agrees well with the Monte Carlo es- 
timates P/iy = 0.04(1) and u ^ 1.2(1) g. In view of 
the above considerations this must be regarded as an ef- 
fective exponent (in our picture the critical transition of 
the gonihedric model lies on the universal ferromagnetic- 
paramagnetic critical surface, and hence has to be de- 
scribed by the usual 3d Ising exponents), induced by a 
crossover phenomenon. Nevertheless, the critical tran- 
sition of the gonihedric model is extremely close to our 
critical end point and this means that the corresponding 
exponents are very good approximations to the critical 
end point ones. Thus from now on we shall use the es- 
timate PcEP = 0.062 ± 0.003. In order to calculate also 
7CEP, we have determined the high temperature suscep- 
tibility for J3 = and J2 = — Ji/4 in the cube approxi- 
mation of the CVM and, according to the CVPAM pre- 
scriptions, we have determined Dlog Fade approximants 
(biased with our improved /3c) to it, and from these we 
have deduced 7cep = 1-41 ± 0.02. 

The analysis of the crossover phenomenon is a con- 
siderably more difficult task, and we have tried to give 
an estimate of the crossover exponent (j) pS] proceeding 
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along the lines described in ||2^. From now on we set 
Jil J\ = Assuming that near the critical end point, 
but still in the ferromagnetic phase, the order parameter 
has a multicritical scaling law given by 

m^^^-P;(,), .^^^^^ (3) 

where t is the deviation from the critical temperature 
and i?cEP can be well approximated by —1/4 (see also 
below), one can derive the scaling laws 

Tc(i?) - Tc(i?cEp) oc (i? - i?cEp)'/'^ (4) 
for the ferromagnetic critical temperature and 

B(i?)oc(i?-i?cEp)-^ ^= /3l-'"S-/^CEP ^ 

where /Rising — 0.327 Q is the usual three-dimensional 
Ising exponent (a high-order CVPAM analysis on the 
simple cubic lattice predicted /3ising = 0.325(4) ||2^), for 
the critical amplitude B{R), which is defined by 

m ~ B{R){T,{R) - T)*--, R > Rcep- (6) 

Using the CVPAM we have extrapolated the low- 
temperature ferromagnetic order parameter (the smallest 
value used being 0.89) given by the CVM star-cube ap- 
proximation, to calculate the critical temperatures and 
amplitudes in the range —0.24 < R < —0.14 for k — 1,2 
and 10. Tc{R) was determined by requiring that Tc{R)- 
biased Dlog Fade approximants gave /3ising = 0.327 and 
then the critical amplitude could be obtained by making 
approximants to {Tc{R) — T)^^^='"sm. As a check, for 
the simple n.n. Ising model (that is k = and i? = 0) 
we have obtained -B(O) ~ 1.626, to be compared with 
the value 1.6919045 reported by Blote and Talapov. The 
results for K = 1 are reported in Tab. |. Several fits 
were then made. The fits on the basis of Eq. ^ gave 
(f) = 1.13,1.09 and 1.03 for k = 1,2 and 10 respectively, 
while the fits to Eq. || gave = 1.16 and 1.20 for k = 1 
and 2 and were inconclusive for k = 10. A fit with i?cEP 
free was also made to Eq. ^ for k = 1 , and the result was 
-Rcep = -0.249774, confirming that i?cEP = -1/4 is a 
very good approximation. 

A reasonable final estimate for the crossover exponent 
might then be = 1.1(1), but it must be taken with 
some care for several reasons. Apart from the various 
approximations involved in the calculations of (j), it is a 
matter of fact that a similar calculation in the case of 
the semi-infinite Ising model |22j gave a result differing 
by 10 to 30% from extensive computer simulations, and 
furthermore the true multicritical scaling law for the or- 
der parameter might be more complicated than our Eq. 
^ (e.g. the scaling axes might not be parallel to the T 
and R axes), although the relatively good quality of the 
fits seems to indicate that Eq. p| is fairly good. 
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Summarizing, we have studied the critical and mul- 
ticritical behavior of the simple cubic Ising model with 
n.n., n.n.n. and plaquette interactions, calculating the 
exponents of a line of critical end points which does not 
exist in two dimensions and might be relevant for some 
surfactant system. We have also tried to give an estimate 
of the crossover exponent, but this can certainly be re- 
fined using e.g. Monte Carlo simulations or series expan- 
sions combined with the partial differential approximants 
method IEtII. 
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TABLE I. Inverse critical temperatures Jc{R) ~ l/Tc{R) and order parameter amplitudes for 
K = 1. 
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FIG. 1. The phase diagram of the model (|l|) for J3 = 0. 
Solid and the dashed lines represent first and second order 
transitions respectively. 



FIG. 2. The phase diagram of the model (^) for k = 1/3. 
Solid and the dashed lines represent first and second order 
transitions respectively. 
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